---
title: "How education, generation and gender jointly structure green and radical right voting"
subtitle: "Reproduction File"
author:
  - name: Armin Schäfer
    affiliation:
      - id: jgu
        name: Johannes Gutenberg University
    orcid: 0000-0002-2538-0092
    email: schaefer@politik.uni-mainz.de
  - name: Nils D. Steiner
    affiliation:
      - ref: jgu
    email: steiner@politik.uni-mainz.de
    orcid: 0000-0003-3433-4079
format: 
  html:
    number-sections: true
    code-fold: true
    code-summary: "Show the code"
    self-contained: true
    toc: true
highlight-style: tango
execute:
  warning: false
  message: false
---

# Preparation

## Packages and appearance

```{r}
#| label: setup
#| echo: true
#| message: false
#| warning: false

if (!require("pacman")) {
  install.packages("pacman")
  library(pacman)
}

p_load(here, haven, sandwich, lme4, tibble, patchwork, pandoc, modelsummary, marginaleffects, gt, tidyverse)

options(
  ggplot2.continuous.colour = "viridis",
  ggplot2.continuous.fill = "viridis",
  dplyr.summarise.inform = FALSE
)

my_theme <- function(...) {
  theme(
    plot.title = element_text(
      color = "black",
      size = 12,
      face = "bold"
    ),
    axis.text = element_text(
      color = "black",
      size = 10
    ),
    axis.title = element_text(
      color = "black",
      face = "bold",
      size = 10
    ),
    panel.grid.minor = element_blank(),
    strip.background = element_rect(
      fill = "white",
      color = "grey70"
    ),
    strip.text = element_text(color = "black"),
    legend.key.size = unit(0.3, "cm")
  )
}
```

## Data and subsets

```{r}
#| label: data
#| code-fold: false

df <- readRDS(here("data", "ess_short.rds"))


dfgreen <- df %>%
  filter(Green_sample == 1 & agea_clean > 24)


dfradicalright <- df %>%
  filter(RR_sample == 1 & agea_clean > 24)

```

# Descriptive results

## Party strength over time

```{r}
#| label: fig-fig_b1
#| fig-cap: "Vote shares of green and radical right parties by country and ESS round (Appendix)"

df %>%
  group_by(cntry, essround) %>%
  summarise(
    mean_rr_vote = mean(vote_radright, na.rm = TRUE),
    mean_green_vote = mean(vote_green, na.rm = TRUE)
  ) %>%
  dplyr::select(essround, mean_rr_vote, mean_green_vote) %>%
  pivot_longer(
    cols = c(mean_rr_vote, mean_green_vote),
    names_to = "party", values_to = "vote"
  ) %>%
  mutate(
    party = case_when(
      party == "mean_rr_vote" ~ "Radical right",
      party == "mean_green_vote" ~ "Green"
    )
  ) %>%
  ggplot(aes(
    x = essround,
    y = vote,
    color = party
  )) +
  geom_line(linewidth = .8) +
  scale_color_manual(values = c("#35b779", "#57106e")) +
  scale_y_continuous(
    labels = scales::percent,
    limits = c(0, .4)
  ) +
  scale_x_continuous(breaks = 1:10) +
  labs(color = NULL, y = "Vote share", x = "Survey round") +
  facet_wrap(~cntry) +
  theme_light() +
  my_theme() +
  theme(
    axis.text = element_text(size = 8),
    strip.text = element_text(size = 8),
    legend.position = "top"
  )

ggsave(
  filename = here("figs", "final", "appendix_fig_b1.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 500
)

```

## Education by generation and gender

```{r}
#| label: fig-fig_b2
#| fig-cap: "Education by cohort and gender (Appendix)"


dfgreen %>%
  filter(!is.na(gndr_fct), !is.na(edu_fct), !is.na(cohort_fct)) %>%
  group_by(cohort_fct, gndr_fct, edu_fct) %>%
  summarise(count = n()) %>%
  mutate(proportion = count / sum(count)) %>%
  ggplot(aes(
    x = fct_relevel(
      cohort_fct,
      "interwar", "post-war",
      "60-70s", "80s", "90s", "millennial"
    ),
    y = proportion,
    fill = fct_relevel(
      edu_fct,
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    )
  )) +
  geom_col(position = "dodge") +
  geom_text(aes(label = scales::percent(proportion, accuracy = 1.0)),
    position = position_dodge(width = 0.9), vjust = -0.5, size = 1.8
  ) +
  labs(
    fill = NULL, x = NULL,
    y = "Proportion of educational degree"
  ) +
  scale_y_continuous(labels = scales::percent) +
  scale_fill_viridis_d(option = "plasma", end = .8) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(
    legend.position = "top",
    axis.text.x = element_text(
      angle = 90,
      vjust = 0.5,
      hjust = 1
    )
  )

ggsave(
  filename = here("figs", "final", "appendix_fig_b2.png"),
  units = c("in"),
  width = 6, height = 4,
  dpi = 300
)

```

## Green vote by generation, education and gender

```{r}
#| label: fig-fig_b3
#| fig-cap: "Voting for a green party by cohort, education, and gender (Appendix)"

dfgreen %>%
  filter(!is.na(gndr_fct), !is.na(edu_fct), !is.na(cohort_fct)) %>%
  group_by(edu_fct, cohort_fct, gndr_fct) %>%
  summarise(mvotegreen = mean(vote_green, na.rm = TRUE)) %>%
  ggplot(aes(
    x = fct_relevel(
      cohort_fct,
      "interwar", "post-war", "60-70s",
      "80s", "90s", "millennial"
    ),
    y = mvotegreen,
    fill = fct_relevel(
      edu_fct,
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    )
  )) +
  geom_col(position = "dodge") +
  geom_text(aes(label = scales::percent(mvotegreen, accuracy = 1.0)),
    position = position_dodge(width = 0.9), vjust = -0.5, size = 2
  ) +
  labs(
    fill = NULL, x = NULL,
    y = "Vote share greens"
  ) +
  scale_y_continuous(labels = scales::percent) +
  scale_fill_viridis_d(option = "plasma", end = .8) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(
    legend.position = "top",
    axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1)
  )

ggsave(
  filename = here("figs", "final", "appendix_fig_b3.png"),
  units = c("in"),
  width = 6, height = 4,
  dpi = 300
)

```

## Radical right vote by generation, education and gender

```{r}
#| label: fig-fig_b4
#| fig-cap: "Voting for a radical right party by cohort, education, and gender (Appendix)"


dfradicalright %>%
  filter(!is.na(gndr_fct), !is.na(edu_fct), !is.na(cohort_fct)) %>%
  group_by(edu_fct, cohort_fct, gndr_fct) %>%
  summarise(mvoterr = mean(vote_radright, na.rm = TRUE)) %>%
  ggplot(aes(
    x = fct_relevel(
      cohort_fct,
      "interwar", "post-war", "60-70s", "80s", "90s", "millennial"
    ),
    y = mvoterr,
    fill = fct_relevel(
      edu_fct,
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    )
  )) +
  geom_col(position = "dodge") +
  geom_text(aes(label = scales::percent(mvoterr, accuracy = 1.0)),
    position = position_dodge(width = 0.9), vjust = -0.5, size = 2
  ) +
  labs(
    fill = NULL, x = NULL,
    y = "Vote share radical right"
  ) +
  scale_y_continuous(labels = scales::percent) +
  scale_fill_viridis_d(option = "plasma", end = .8) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(
    legend.position = "top",
    axis.text.x = element_text(angle = 90, vjust = 0.5, hjust = 1)
  )

ggsave(
  filename = here("figs", "final", "appendix_fig_b4.png"),
  units = c("in"),
  width = 6, height = 4,
  dpi = 300
)

```

# Green vote

## Multilevel model: Voting for the greens (two-way interaction)

```{r}
#| label: fig-fig_1
#| fig-cap: "Voting for a green party by cohort and education (twoway interaction)"
#| fig-width: 8
#| fig-height: 8

# model
two_green <- glmer(
  vote_green ~
    cohort_fct * edu_fct +
    age_grp_fct * edu_fct +
    factor(essround) * edu_fct +
    gndr_fct + immigrant + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfgreen
)

# predicted probabilities
pp_two_green <- predictions(two_green,
  by = c("cohort_fct", "edu_fct"),
  newdata = datagrid(
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
    ),
    re.form = NA
)

plot_pp_two_green <- ggplot(
  data = pp_two_green,
  aes(
    x = cohort_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL,
    y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_1_upper.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# marginal effects
me_two_green <- avg_slopes(
  two_green,
  variables = "edu_fct",
  by = c("cohort_fct"),
  newdata = datagrid(
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
  re.form=NA
)

me_two_green <- me_two_green %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_two_green <- ggplot(data = me_two_green, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = contrast
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_1_lower.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# Combine the plots vertically
combined_plot <- plot_pp_two_green / plot_me_two_green

# Save the combined plot
ggsave(
  filename = here("figs", "final", "figure_1.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10,
  dpi = 300
)

combined_plot

```

## Voting for a green party -- education and ESS wave

```{r}
#| label: fig-fig_d1
#| fig-cap: "Predicted probabilities of voting for a green party by education and ESS wave (Appendix)"
#| fig-width: 8

# predicted probabilities: by period
pp_two_green_period <- predictions(
  two_green,
  by = c("essround", "edu_fct"),
  newdata = datagrid(
    essround = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_two_green_period,
  aes(
    x = essround, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d1.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Voting for a green party -- education and age group

```{r}
#| label: fig-fig_d2
#| fig-cap: "Predicted probabilities of voting for a green party by education and age group (Appendix)"
#| fig-width: 8

# predicted probabilities: by age group
pp_two_green_age <- predictions(
  two_green,
  by = c("age_grp_fct", "edu_fct"),
  newdata = datagrid(
    age_grp_fct = c("22-29", "30-65", "above 65"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_two_green_age,
  aes(
    x = age_grp_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d2.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)
```

## Multilevel model: Voting for the greens (three-way interaction)

```{r}
#| label: fig-fig_2
#| fig-cap: "Predicted probabilities and marginal effects of voting for a green party (three-way interaction)"
#| fig-width: 8
#| fig-height: 8


# model
three_green <- glmer(
  vote_green ~
    cohort_fct * edu_fct * gndr_fct +
    age_grp_fct * edu_fct * gndr_fct +
    factor(essround) * edu_fct * gndr_fct +
    immigrant + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfgreen
)

# predicted probabilities
pp_three_green <- predictions(three_green,
  by = c("gndr_fct", "cohort_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

plot_pp_three_green <- ggplot(
  data = pp_three_green,
  aes(
    x = cohort_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "fig_2_upper.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# marginal effects
me_three_green_edu <- avg_slopes(
  three_green,
  variables = "edu_fct",
  by = c("cohort_fct", "gndr_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_green_edu <- me_three_green_edu %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_three_green <- ggplot(data = me_three_green_edu, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = contrast
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_2_lower.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# Combine the plots vertically
combined_plot <- plot_pp_three_green / plot_me_three_green

# Save the combined plot
ggsave(
  filename = here("figs", "final", "figure_2.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10,
  dpi = 300
)

combined_plot

```

## Table with predicted probabilities of voting for the greens

```{r}
#| label: tbl-table_pp_green
#| tbl-cap: "Predicted probabilities of voting for the greens (threeway interaction)"

pp_three_green %>% 
  select(cohort_fct:p.value) %>% 
  rename(cohort = cohort_fct) %>% 
gt() %>% 
  fmt_number(
    decimals = 3
  )

```

## Multilevel model: Voting for the greens -- education, ESS round & gender (three-way interaction)

```{r}
#| label: fig-fig_d3
#| fig-cap: "Predicted probabilities of voting for a green party by education, ESS wave and gender (Appendix)"
#| fig-width: 8


# predicted probabilities: by period
pp_three_green_period <- predictions(three_green,
  by = c("gndr_fct", "essround", "edu_fct"),
  newdata = datagrid(
    essround = c("1", "2", "3", "4", "5", "6", "7", "8", "9", "10"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_three_green_period,
  aes(
    x = essround,
    y = estimate,
    ymin = conf.low,
    ymax = conf.high,
    color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "appendix_fig_d3.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Multilevel model: Voting for the greens -- education, ESS round & age group (three-way interaction)

```{r}
#| label: fig-fig_d4
#| fig-cap: "Predicted probabilities of voting for a green party by education, ESS wave and age group (Appendix)"
#| fig-width: 8


# predicted probabilities: by age_group
pp_three_green_age <- predictions(three_green,
  by = c("gndr_fct", "age_grp_fct", "edu_fct"),
  newdata = datagrid(
    age_grp_fct = c("22-29", "30-65", "above 65"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_three_green_age,
  aes(
    x = age_grp_fct,
    y = estimate,
    ymin = conf.low,
    ymax = conf.high,
    color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "appendix_fig_d4.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Marginal effects of women vs. men on voting for the greens

```{r}
#| label: fig-fig_d5
#| fig-cap: "Marginal effects of women vs. men on voting for the greens across cohorts and education levels (Appendix)"
#| fig-width: 8

# contrasting women and men
me_three_green_gender <- avg_slopes(
  three_green,
  variables = "gndr_fct",
  by = c("cohort_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_green_gender <- me_three_green_gender %>%
  mutate(contrast = case_when(
    contrast == "women - men" ~
      "women vs men"
  ))

ggplot(data = me_three_green_gender, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = edu_fct
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of women vs. men"
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d5.png"),
  units = c("in"),
  width = 6, height = 5,
  dpi = 300
)
```

## Remove files

```{r}
#| label: remove_files

rm(list = c("df", "three_green", "two_green", "me_three_green_edu", "me_three_green_gender",
            "me_two_green", "pp_two_green", "pp_two_green_age", "pp_two_green_period",
               "pp_three_green", "pp_three_green_age", "pp_three_green_period"))

```

# Radical right vote

## Multilevel model: Voting for the radical right (two-way interaction)

```{r}
#| label: fig-fig_3
#| fig-cap: "Voting for a radical right party by cohort and education (twoway interaction)"
#| fig-width: 8
#| fig-height: 8

# model
two_rr <- glmer(
  vote_radright ~
    cohort_fct * edu_fct +
    age_grp_fct * edu_fct +
    factor(essround) * edu_fct +
    gndr_fct + immigrant + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfradicalright
)

# predicted probabilities
pp_two_rr <- predictions(two_rr,
  by = c("cohort_fct", "edu_fct"),
  newdata = datagrid(
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

plot_pp_two_rr <- ggplot(
  data = pp_two_rr,
  aes(
    x = cohort_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_3_upper.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# marginal effects
me_two_rr <- avg_slopes(
  two_rr,
  variables = "edu_fct",
  by = c("cohort_fct"),
  newdata = datagrid(
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_two_rr <- me_two_rr %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_two_rr <- ggplot(data = me_two_rr, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = contrast
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_3_lower.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# Combine the plots vertically
combined_plot <- plot_pp_two_rr / plot_me_two_rr

# Save the combined plot
ggsave(
  filename = here("figs", "final", "figure_3.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10,
  dpi = 300
)

combined_plot

```

## Voting for a radical right party by education and ESS round

```{r}
#| label: fig-fig_d6
#| fig-cap: "Predicted probabilities of voting for a radical right party by education and ESS wave (Appendix)"
#| fig-width: 8

# predicted probabilities: by period
pp_two_rr_period <- predictions(two_rr,
  by = c("essround", "edu_fct"),
  newdata = datagrid(
    essround = c(
      "1",
      "2",
      "3",
      "4",
      "5",
      "6",
      "7",
      "8",
      "9",
      "10"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_two_rr_period,
  aes(
    x = essround, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d6.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Voting for a radical right party by education and age group

```{r}
#| label: fig-fig_d7
#| fig-cap: "Predicted probabilities of voting for a radical right party by education and age group (Appendix)"
#| fig-width: 8

# predicted probabilities: by age
pp_two_rr_age <- predictions(two_rr,
  by = c("age_grp_fct", "edu_fct"),
  newdata = datagrid(
    age_grp_fct = c("22-29", "30-65", "above 65"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_two_rr_age,
  aes(
    x = age_grp_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d7.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)
```

## Multilevel model: Voting for radical right parties (three-way interaction)

```{r}
#| label: fig-fig_4
#| fig-cap: "Predicted probabilities and marginal effects of voting for a radical right party"
#| fig-width: 8
#| fig-height: 8


# model
three_rr <- glmer(
  vote_radright ~ cohort_fct * edu_fct * gndr_fct + age_grp_fct * edu_fct * gndr_fct +
    immigrant + factor(essround) * edu_fct * gndr_fct + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfradicalright
)

# predicted probabilities
pp_three_rr <- predictions(three_rr,
  by = c("gndr_fct", "cohort_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

plot_pp_three_rr <- ggplot(
  data = pp_three_rr,
  aes(
    x = cohort_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "fig_4_upper.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# marginal effects
me_three_rr_edu <- avg_slopes(
  three_rr,
  variables = "edu_fct",
  by = c("cohort_fct", "gndr_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_rr_edu <- me_three_rr_edu %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_three_rr <- ggplot(data = me_three_rr_edu, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = contrast
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "fig_4_lower.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

# Combine the plots vertically
combined_plot <- plot_pp_three_rr / plot_me_three_rr

# Save the combined plot
ggsave(
  filename = here("figs", "final", "figure_4.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10,
  dpi = 300
)

combined_plot

```

## Table with predicted probabilities of voting for the radical right

```{r}
#| label: tbl-table_pp_rr
#| tbl-cap: "Predicted probabilities of voting for the radical right (threeway interaction)"

pp_three_rr %>% 
  select(cohort_fct:p.value) %>% 
  rename(cohort = cohort_fct) %>% 
gt() %>% 
  fmt_number(
    decimals = 3
  )

```

## Predicted probabilities of radical right: by education, ESS round and gender

```{r}
#| label: fig-fig_d8
#| fig-cap: "Predicted probabilities of voting for a radical right party by education, ESS wave and gender (Appendix)"
#| fig-width: 8

pp_three_rr_period <- predictions(three_rr,
  by = c("gndr_fct", "essround", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    essround = c(
      "1",
      "2",
      "3",
      "4",
      "5",
      "6",
      "7",
      "8",
      "9",
      "10"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_three_rr_period,
  aes(
    x = essround, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "appendix_fig_d8.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Predicted probabilities of radical right: by gender, cohort and age group

```{r}
#| label: fig-fig_d9
#| fig-cap: "Predicted probabilities of voting for a radical right party by education, age group and gender (Appendix)"
#| fig-width: 8


pp_three_rr_age <- predictions(three_rr,
  by = c("gndr_fct", "age_grp_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    age_grp_fct = c("22-29", "30-65", "above 65"),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

ggplot(
  data = pp_three_rr_age,
  aes(
    x = age_grp_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

ggsave(
  filename = here("figs", "final", "appendix_fig_d9.png"),
  units = c("in"),
  width = 8, height = 5,
  dpi = 300
)

```

## Marginal effects of women vs. men on voting for the radical right across cohorts and education levels

```{r}
#| label: fig-fig_d10
#| fig-cap: "Marginal effects of women vs. men on voting for the radical right across cohorts and education levels (Appendix)"
#| fig-width: 8

me_three_rr_gender <- avg_slopes(
  three_rr,
  variables = "gndr_fct",
  by = c("cohort_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_fct = c(
      "interwar",
      "post-war",
      "60-70s",
      "80s",
      "90s",
      "millennial"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_rr_gender <- me_three_rr_gender %>%
  mutate(contrast = case_when(
    contrast == "women - men" ~
      "women vs men"
  ))

ggplot(data = me_three_rr_gender, aes(
  x = cohort_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = edu_fct
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of women vs. men"
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

ggsave(
  filename = here("figs", "final", "appendix_fig_d10.png"),
  units = c("in"),
  width = 6, height = 5,
  dpi = 300
)
```

# Tables

## Preliminaries

```{r}
#| label: table_setup

coefname <- c(
  "cohort_fctpost-war" = "Post-War Generation",
  "cohort_fct60-70s" = "60-70s Generation",
  "cohort_fct80s" = "80s Generation",
  "cohort_fct90s" = "90s Generation",
  "cohort_fctmillennial" = "Millennials",
  "gndr_fctwomen" = "Gender: female",
  "factor(age_group)2" = "Age: 30-65",
  "factor(age_group)3" = "Age: 66+",
  "immigrant" = "Immigrant background",
  "edu_fctupper secondary" = "Education: upper secondary",
  "edu_fcttertiary" = "Education: tertiary",

  # cohort # education
  "cohort_fctpost-war:edu_fctupper secondary" =
    "Post-War Generation # education: upper secondary",
  "cohort_fctpost-war:edu_fcttertiary" =
    "Post-War Generation # education: tertiary",
  "cohort_fct60-70s:edu_fctupper secondary" =
    "60-70s Generation # education: upper secondary",
  "cohort_fct60-70s:edu_fcttertiary" =
    "60-70s Generation # education: tertiary",
  "cohort_fct80s:edu_fctupper secondary" =
    "80s Generation # education: upper secondary",
  "cohort_fct80s:edu_fcttertiary" =
    "80s Generation # education: tertiary",
  "cohort_fctmillennial:edu_fctupper secondary" =
    "Millennials # education: upper secondary",
  "cohort_fctmillennial:edu_fcttertiary" =
    "Millennials # education: tertiary",

  # cohort # gender
  "cohort_fctpost-war:gndr_fctwomen" =
    "Post-War Generation # Female",
  "cohort_fct60-70s:gndr_fctwomen" =
    "60-70s Generation # Female",
  "cohort_fct80s:gndr_fctwomen" =
    "80s Generation # Female",
  "cohort_fct90s:gndr_fctwomen" =
    "90s Generation # Female",
  "cohort_fctmillennial:gndr_fctwomen" =
    "Millennials # Female",

  # education # gender
  "edu_fctupper secondary:gndr_fctwomen" =
    "Education: upper secondary # Female",
  "edu_fcttertiary:gndr_fctwomen" =
    "Education: tertiary # female",

  # cohort # education # gender
  "cohort_fctpost-war:edu_fctupper secondary:gndr_fctwomen" =
    "Post-War Generation # Educ.: upper secondary # Female",
  "cohort_fct60-70s:edu_fctupper secondary:gndr_fctwomen" =
    "60-70s Generation # Educ.: upper secondary # Female",
  "cohort_fct80s:edu_fctupper secondary:gndr_fctwomen" =
    "80s Generation # Educ.: upper secondary # Female",
  "cohort_fct90s:edu_fctupper secondary:gndr_fctwomen" =
    "90s Generation # Educ.: upper secondary # Female",
  "cohort_fctmillennial:edu_fctupper secondary:gndr_fctwomen" =
    "Millennials # Educ.: upper secondary # Female",
  "cohort_fctpost-war:edu_fcttertiary:gndr_fctwomen" =
    "Post-War Generation # Educ.: tertiary # Female",
  "cohort_fct60-70s:edu_fcttertiary:gndr_fctwomen" =
    "60-70s Generation # Educ.: tertiary # Female",
  "cohort_fct80s:edu_fcttertiary:gndr_fctwomen" =
    "80s Generation # Educ.: tertiary # Female",
  "cohort_fct90s:edu_fcttertiary:gndr_fctwomen" =
    "90s Generation # Educ.: tertiary # Female",
  "cohort_fctmillennial:edu_fcttertiary:gndr_fctwomen" =
    "Millennials # Educ.: tertiary # Female"
)


# add rows
rows <- tribble(
  ~term, ~"Model 1", ~"Model 2", ~"Model 3", ~"Model 4",
  "Interwar Generation", "ref.", "ref.", "ref.", "ref.",
  "Gender: male", "ref.", "ref.", "ref.", "ref.",
  "Age: 22-29", "ref.", "ref.", "ref.", "ref.",
  "Education: (Below) lower secondary", " ", "ref.", "ref.", "ref.",
  "Generation X Education", "", "", "", "",
  "Generation X Female", "", "", "", "",
  "Education X Female", "", "", "", "",
  "Generation X Education X Female", "", "", "", "",
  "Country fixed effects", "yes", "yes", "yes", "yes",
  "Twoway interaction Age group X Education", "no", "no", "yes", "yes",
  "Twoway interaction ESS wave X Education", "no", "no", "yes", "yes",
  "Threeway interaction Age group X Education X Gender", "no", "no", "no", "yes",
  "Threeway interaction ESS wave X Education X Gender", "no", "no", "no", "yes"
)

# position of additional rows
attr(rows, "position") <- c(1, 12, 15, 22, 27, 44, 55, 60, 81, 82, 83, 84, 85)


gof <- c("nobs", "r.squared", "AIC", "BIC")
```

## All models of green vote

```{r}
#| label: tbl-table_b1
#| fig-cap: "Voting for green parties -- regression coefficients"


models <- list(
  "Model 1" = glmer(
    vote_green ~ cohort_fct + gndr_fct +
      +factor(age_group) +
      factor(essround) + factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfgreen
  ),
  "Model 2" = glmer(
    vote_green ~ cohort_fct + gndr_fct + edu_fct +
      factor(age_group) + immigrant +
      factor(essround) + factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfgreen
  ),
  "Model 3" = glmer(
    vote_green ~ cohort_fct * edu_fct +
      factor(age_group) * edu_fct +
      factor(essround) * edu_fct +
      gndr_fct + immigrant +
      +factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfgreen
  ),
  "Model 4" = glmer(
    vote_green ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct +
      immigrant +
      factor(cntry_id) + (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfgreen
  )
)


# create the table
modelsummary(models,
  stars = TRUE, fmt = 2,
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = here("tables", "appendix_table_c1.docx")
)

modelsummary(models,
  stars = TRUE, fmt = 2,
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = "gt")

```

## All models of radical right vote

```{r}
#| label: tbl-table_b2
#| fig-cap: "Voting for a radical right parties -- regression coefficients"


models <- list(
  "Model 1" = glmer(
    vote_radright ~ cohort_fct + gndr_fct +
      factor(age_group) +
      factor(essround) + factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfradicalright
  ),
  "Model 2" = glmer(
    vote_radright ~ cohort_fct + gndr_fct + edu_fct +
      factor(age_group) + immigrant +
      factor(essround) + factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfradicalright
  ),
  "Model 3" = glmer(
    vote_radright ~ cohort_fct * edu_fct +
      factor(age_group) * edu_fct +
      factor(essround) * edu_fct +
      gndr_fct + immigrant +
      +factor(cntry_id) +
      (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfradicalright
  ),
  "Model 4" = glmer(
    vote_radright ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct +
      immigrant +
      factor(cntry_id) + (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfradicalright
  )
)


# create the table
modelsummary(models,
  stars = TRUE, fmt = 2,
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = here("tables", "appendix_table_c2.docx")
)

modelsummary(models,
  stars = TRUE, fmt = 2,
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = "gt")

```

# Alternative models

## Logistic regression with country-wave fixed effects and clustered standard errors.

```{r}
#| label: tbl-table_b3
#| fig-cap: "Voting for a radical right parties and green parties--multilevel and logitistic regression models in comparison"

models <- list(
  "Multilevel Greens" = glmer(
    vote_green ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct +
      immigrant +
      factor(cntry_id) + (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfgreen
  ),
  "Logit Greens" = glm(
    vote_green ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct + immigrant +
      factor(cntry_wave_id),
    family = binomial(link = "logit"), data = dfgreen
  ),
  "Multilevel Radical Right" = glmer(
    vote_radright ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct +
      immigrant + factor(cntry_id) + (1 | cntry_wave_id),
    family = binomial(link = "logit"),
    nAGQ = 0,
    data = dfradicalright
  ),
  "Logit Radical Right" = glm(
    vote_radright ~ cohort_fct * edu_fct * gndr_fct +
      factor(age_group) * edu_fct * gndr_fct +
      factor(essround) * edu_fct * gndr_fct + immigrant +
      factor(cntry_wave_id),
    family = binomial(link = "logit"), data = dfradicalright
  )
)

rows <- tribble(
  ~term, ~"Model 1", ~"Model 2", ~"Model 3", ~"Model 4",
  "Interwar Generation", "ref.", "ref.", "ref.", "ref.",
  "Gender: male", "ref.", "ref.", "ref.", "ref.",
  "Age: 22-29", "ref.", "ref.", "ref.", "ref.",
  "Education: (Below) lower secondary", " ", "ref.", "ref.", "ref.",
  "Generation X Education", "", "", "", "",
  "Generation X Female", "", "", "", "",
  "Education X Female", "", "", "", "",
  "Generation X Education X Female", "", "", "", "",
  "Country fixed effects", "yes", "", "yes", "",
  "Country-wave fixed effects", "", "yes", "", "yes",
  "Threeway interaction Age group X Education X Gender", "yes", "yes", "yes", "yes",
  "Threeway interaction ESS wave X Education X Gender", "yes", "yes", "yes", "yes"
)

attr(rows, "position") <- c(1, 12, 15, 22, 27, 44, 55, 60, 81, 82, 83, 84)

modelsummary(models,
  stars = TRUE, fmt = "%.3f",
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = here("tables", "appendix_table_e1.docx")
)

modelsummary(models,
  stars = TRUE, fmt = "%.3f",
  coef_map = coefname,
  add_rows = rows,
  gof_omit = "ICC|RMSE",
  output = "gt")

```

# Fine-grained cohorts

## Multilevel model with three-way interaction: Voting for the greens (alternative cohorts)

```{r}
#| label: fig-fig_f1
#| fig-cap: "Voting for a green party by cohort and education (three-way interaction)"
#| fig-width: 8
#| fig-height: 8

# model
three_green_alt <- glmer(
  vote_green ~
    cohort_alt_fct * edu_fct * gndr_fct +
    factor(age_group) * edu_fct * gndr_fct +
    factor(essround) * edu_fct * gndr_fct +
    immigrant + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfgreen
)

# predicted probabilities
pp_three_green_alt <- predictions(three_green_alt,
  by = c("gndr_fct", "cohort_alt_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_alt_fct = c(
      "-1929",
      "1930s",
      "1940s",
      "1950s",
      "1960s",
      "1970s",
      "1980s",
      "1990+"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

# plot
plot_pp_three_green_alt <- ggplot(
  data = pp_three_green_alt,
  aes(
    x = cohort_alt_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

# marginal effects
me_three_green_alt <- avg_slopes(
  three_green_alt,
  variables = "edu_fct",
  by = c("cohort_alt_fct", "gndr_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_alt_fct = c(
      "-1929",
      "1930s",
      "1940s",
      "1950s",
      "1960s",
      "1970s",
      "1980s",
      "1990+"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_green_alt <- me_three_green_alt %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_three_green_alt <- ggplot(
  data = me_three_green_alt,
  aes(
    x = cohort_alt_fct,
    y = estimate,
    ymin = conf.low,
    ymax = conf.high,
    color = contrast
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL,
    color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

# Combine the plots vertically
combined_plot <- plot_pp_three_green_alt / plot_me_three_green_alt

# Save the combined plot
ggsave(
  filename = here("figs", "final", "appendix_fig_f1.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10, # Adjust height to fit two plots
  dpi = 300
)

combined_plot

```

## Multilevel model with three-way interaction: Voting for the radical right (alternative cohorts)

```{r}
#| label: fig-fig_f2
#| fig-cap: "Voting for a green party by cohort and education (three-way interaction)"
#| fig-width: 8
#| fig-height: 8

# model
three_rr_alt <- glmer(
  vote_radright ~
    cohort_alt_fct * edu_fct * gndr_fct +
    factor(age_group) * edu_fct * gndr_fct +
    factor(essround) * edu_fct * gndr_fct +
    immigrant + factor(cntry_id) +
    (1 | cntry_wave_id),
  family = binomial(link = "logit"),
  nAGQ = 0,
  data = dfradicalright
)

# predicted probabilities
pp_three_rr_alt <- predictions(three_rr_alt,
  by = c("gndr_fct", "cohort_alt_fct", "edu_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_alt_fct = c(
      "-1929",
      "1930s",
      "1940s",
      "1950s",
      "1960s",
      "1970s",
      "1980s",
      "1990+"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)


# plot
plot_pp_three_rr_alt <- ggplot(
  data = pp_three_rr_alt,
  aes(
    x = cohort_alt_fct, y = estimate, ymin = conf.low,
    ymax = conf.high, color = edu_fct
  )
) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, y = "Predicted probability",
    color = NULL
  ) +
  scale_color_viridis_d(option = "plasma", end = .8) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top") +
  facet_wrap(~gndr_fct)

# marginal effects
me_three_rr_alt <- avg_slopes(
  three_rr_alt,
  variables = "edu_fct",
  by = c("cohort_alt_fct", "gndr_fct"),
  newdata = datagrid(
    gndr_fct = c("men", "women"),
    cohort_alt_fct = c(
      "-1929",
      "1930s",
      "1940s",
      "1950s",
      "1960s",
      "1970s",
      "1980s",
      "1990+"
    ),
    edu_fct = c(
      "(below) lower secondary",
      "upper secondary",
      "tertiary"
    ),
    grid_type = "counterfactual"
  ),
    re.form = NA
)

me_three_rr_alt <- me_three_rr_alt %>%
  mutate(contrast = case_when(
    contrast == "tertiary - (below) lower secondary" ~
      "high vs low edu",
    contrast == "upper secondary - (below) lower secondary" ~
      "medium vs low edu"
  ))

plot_me_three_rr_alt <- ggplot(data = me_three_rr_alt, aes(
  x = cohort_alt_fct, y = estimate,
  ymin = conf.low, ymax = conf.high,
  color = contrast
)) +
  geom_hline(yintercept = 0, color = "grey70") +
  geom_pointrange(position = position_dodge(width = c(0.2))) +
  labs(
    x = NULL, color = "Education",
    y = "Marginal effect of education"
  ) +
  scale_color_manual(values = c("#fe9f6d", "#8c2981")) +
  facet_wrap(~gndr_fct) +
  theme_light() +
  my_theme() +
  theme(legend.position = "top")

# Combine the plots vertically
combined_plot <- plot_pp_three_rr_alt / plot_me_three_rr_alt

# Save the combined plot
ggsave(
  filename = here("figs", "final", "appendix_fig_f2.png"),
  plot = combined_plot,
  units = c("in"),
  width = 8,
  height = 10, # Adjust height to fit two plots
  dpi = 300
)

combined_plot

```
